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Abstract 


In  this  paper  we  show  that  the  M/M/c  queue,  with 
arrival  and  service  rates  which  vary  according  to  the 
state  of  a Markov  process,  has  a steady-state  probability 
vector  of  a modified  matrix-geometric  form.  The  rate 
matrix  R is  the  unique  positive  solution  to  a quadratic 
matrix  equation,  which  may  be  solved  numerically  by  suc- 
cessive substitutions.  A theorem  which  provides  an  accur- 
acy check  on  that  computation  is  proved. 

Finally  a numerical  example  is  discussed  and  its 
results  are  interpreted. 


. Key  Words 

Queueing  theory,  M/M/c  queue,  random  environment, 
fluctuating  queues,  computational  probability. 
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I.  Introduction 


This  paper  contains  generalizations  and  further 
elaborations  of  results  discussed  in  [1] . We  refer  to  this 
earlier  paper  for  all  definitions  and  notations.  The  main 
new  topics  of  discussion  are:  the  joint  stationary  proba- 
bility distribution  of  the  queue  length  and  the  underlying 
phase  state,  considered  at  epochs  where  the  phase  state 
changes;  the  extension  of  Theorem  5 of  [1]  to  the  M/M/c 
queue;  the  discussion  of  some  computational  aspects  and 
numerical  examples.  Throughout  this  paper,  we  shall  only 
consider  stable  versions  of  the  queues  under  discussion. 

Lemma  1 

Provided  that  tt_X  < ru,  the  positive  matrix  R,  which 
is  the  unique  solution  to  the  equation 

(1)  R2My)  + R[Q-A(A)-Mji>  ] + MX)  =0, 

in  the  set  of  nonnegative  matrices  of  spectral  radius  less 
than  one,  satisfies 

(2)  R u = X. 

Proof 

The  matrix  R was  introduced  in  Theorem  5 of  [1] • By 
postmultiplying  (1)  by  e,  we  obtain 
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R2£  - RU+£)  + X = ( I-R)  U-Rp)  = 0. 

Since  the  spectral  radius  of  R is  less  than  one,  the  result 
follows . 

Remarks 

a.  The  equality  (2)  can  either  be  used  as  an  accuracy  check 
on  the  numerical  computation  of  R or  can  be  incorporated  into 
the  evaluation  of  R,  so  as  to  expedite  that  algorithm. 

b.  The  inner  product  irRe  is  easily  seen  to  be  the  steady- 
state  probability  that  the  server  is  occupied.  It  should  be 
noted  that  in  general  this  quantity  is  not  equal  to 

p = (irp) -1.  In  pedagogical  material,  it  should  therefore 

be  stressed  that  the  interpretation  given  to  the  traffic 
intensity  is  strongly  model-dependent. 

An  exception  is  the  interesting  special  case,  where 
£ = p*e.  In  that  case,  the  equality  p = jiRe  follows  readily 
from  (2) . 
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II.  The  Queue  Length  at  Phase  Transitions 

The  motivation  for  considering  this  embedded  process 
is  of  some  interest.  As  noted  in  [ 1 ] , the  most  interesting 
cases  of  the  M/M/1  queue  in  a Markovian  environment  involve 
one  or  more  phases  j during  which  Aj  > yj,  even  though  the 
global  queue  is  stable.  The  realizations  of  such  queueing 
processes  exhibit  substantial  random  oscillations  which  are 
evident  in  their  effect  on  the  components  of  the  steady- 
state  vector  x of  Thm.  5 [ 1 ] . During  a sojourn  in  a phase 
j for  which  the  queue  is  "locally"  unstable,  substantial 
buildup  may  occur  to  be  cleared  during  later  phases  for  which 
the  queue  is  locally  stable. 

One  might  expect  to  obtain  useful  additional  informa- 
tion by  considering  the  queue  lengths  at  the  successive  phase 
transitions.  We  shall  derive  the  joint  stationary  probabil- 
ity density  of  the  queue  length  and  the  phase  immediately 
before  and  immediately  after  phase  transitions.  We  shall 
however  show  that  the  conditional  queue  length  density  at  the 
end  of  a sojourn  in  the  phase  j is  the  same  as  the  conditional 
queue  length  density,  given  that  the  phase  is  j . By  suitable 
interpretation  of  a detailed  numerical  example  in  Section  IV, 
we  shall  see  that  this  result  is  due  to  the  exponential 
nature  of  the  sojourn  times  in  the  various  phases. 
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It  is  clear  that  the  sequence  { ( 5n, Jn) ,n>0 } , where 
5n  is  the  queue  length  and  Jn  the  phase  state  .immediately 
after  the  n-th  phase  transition,  is  a Markov  chain  on  the 
state  space  { (i,h) : ilO,  lfhim).  We  denote  by  p^^i(h,u), 
the  probability  that  in  an  ordinary  M/M/1  queue  with  para- 
meters A^,  the  queue  length  at  time  u is  i',  given  that 

it  is  i at  time  0.  Explicit  formulas  for  pj^,(h,u)  are 
known,  but  will  not  be  needed  in  our  discussion.  The  transi- 
tion probability  matrix  P of  the  chain  HSn,Jn)}  is  then 
given  by 

g,  _ahU 

(3)  Phh,(i,i,)  = 0^°°Pii'  (h'u)e  duQhh’' 

for  iiO,  i'iO,  h^h' , lih,  h'<m. 

= 0,  for  i>0,  i' >0,  h=h' . 

The  exponent  ah  = -Qhh,  for  1-h-m.  We  shall  find  it 
convenient  to  partition  the  matrix  P into  mxm  blocks 
P(i,i')  = {Phh1 (i • i ' ) J , which  we  write  as 

(4)  P ( i , i ' ) = V(i,i' ) [Q+A (a) ] , for  i>0,  i'>0, 

where  A (a)  = diag (a^, . . . ,am)  and  V(i,i') 

.0° 

= diag  / p. . , (h,u) e 

lo 

Theorem  1 

The  stationary  probability  vector  £ = (y.0'Xl'’**) 
of  P is  given  by 


“°hu. 


Lfhimj 
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(5)  Y.i  ~ (l£)  ^(I-R)  [Q+A  (a_)  ] , for  i>0. 

Proof 

The  steady-state  equations  y = yP,  may  be  written  as 

oo 

(6)  Zi  = [ L.V(i'fi)[Q+4(o)]  =wi[Q+A(a)], 

i ' =0 


for  i>0,  where  w^  = 7 y. ,V(i',i). 

i'=0 


Setting  £ v.  = u,  we  first  obtain  that 
i=0 


-1 

uA  (£) [Q+A(£)]=u.  Also  ue=l , since  ye=l.  It  readily  follows 
that  u = (ttct)  A (ct)  , a result  which  is  to  be  expected. 


-o,  u 


— Ul  u 

Setting  / p.,.(h,u)e  du  = p*  . (h) , for  l<h<m,  the 
0 1 1 11 

birth-and-death  equations  for  the  M/M/1  queue  yield  in  a 
straightforward  manner  that 


(?)  pi*o (h) (ah+Ah)  * 6i'0  + Pi'l<h>Wh 

pi«i(h)  (ah+xh+V  = 6i'i  + pi',i-l(h)xh  + pi'  ,i+i (h)  uh' 

for  i'10,  i>l,  l<him. 

It  now  readily  follows  that 
(8)  V(i*  ,0)  a(£+X)  - $it()I  + V(i\l)  A(u)  , 

V(i’  ,i)  A (o+X+£)  = 6±  , ±I  + V(i' ,i-l) A(X) 

+ V(i' ,i+l) A(y)  , 

for  i'20,  i>l. 
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We  further  obtain  that 

w0A(a  + X)  = y.0+w1A(ja)  = Wg  [Q+A  (a)  ] +w^A  (_^) 

w^A  (o+X+u)  = Jli+wi_1A  ( X ) +wi+1A  (y_) 

= w±  [Q+A  (aj  ] +W  ^ A (A.)  +Wj_+^A  (_y)  , 

for  i>0. 

Upon  simplification,  we  notice  that  the  vectors 
wn-  , ilO , satisfy  the  equations 

(10)  Wg[Q-A(X)]  + W1A(y)  = 0 

-i-lA^-)  + ^itQ-A (A+£)  1 + wi+1A(u)  = £, 

for  ill.  These  are,  except  for  the  normalizing  condition, 

precisely  the  same  equations  as  satisfied  by  the  vectors 

x^,  ilO,  in  Thm.  5 of  [1],  Moreover,  the  vectors  need  to 

be  positive  and  the  sum  ? w • e must  be  finite.  It  now  fol- 

i=0 

lows  from  classical  results  an  irreducible,  positive  recur- 
rent Markov  chains  that  w^  = kx^,  for  ilO  and  some  positive 
constant  k.  It  follows  immediately  that  ^=k^( I-R)  R1  [Q+A  ( a^)  ] , 
for  ilO,  which  upon  normalization  yields  the  stated  result. 

Corollary  1 

The  conditional  mean  queue  length  at  the  beginning 
of  a j -phase  is  given  by 

(11)  — 1 ..Jit(I-R)  -1R  [Q+A  ( a ) ]l  , for  l<j<m. 

1Tjaj[  fj 

The  conditional  stationary  density  of  the  queue 
length  at  the  beginning  of  a j -phase  is  given  by 
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(12)  — i — i 7T  ( I-R)  R1  [ Q+A  (a ) ] i , for  i>0. 

1 ~ " ]j 

Corollary  2 

Let  z = (Zq,^,...)  with  components  z^j  be  the 
steady-state  probability  vector  of  the  queue  length  and 
the  phase  immediately  before  a phase  change,  then 

(13)  z^  = (U£)  ■'"ir  (I-R)  R^A  {a)  , for  ilO. 

Proof 

The  vectors  z^  and  ^ are  related  by 

(14)  Zj_A  1(a)[Q+A(a)]  = , for  iiO. 

Clearly  the  vectors  z. in  (13)  are  positive,  satisfy  the 
equations  (14)  and  ze-1.  This  guarantees  that  z_  is  the 
steady-state  vector  of  the  Markov  chain,  obtained  by  con- 
sidering the  queue  length  and  phase  immediately  prior  to 
phase  transitions. 

Remark 

We  see  that 

(15)  1Ijlxij  = (Z2.)  (tfjOj)  for  i-0'  l-j-m, 

so  that  the  conditional  probabilities  that  there  are  i 
customers  in  the  queue,  given  that  the  phase  is  j and  given 
that  a sojourn  in  the  phase  j has  ended,  are  equal. 
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III.  The  Multi-server  Queue 

The  M/M/c  queue  with  randomly  varying  arrival  and 
service  rates  defines  a Markov  process  with  the  infinitesimal 
generator  Q*(c),  given  by 


(15)  Q*(c) 


Aoo 

A01 

A10 

A11 

A12 

A20 

A21 

CM 

CM 

< 

Ac-1 , 0 Ac-l,l  Ac-1, 2 

A0  A2 

An  A. 


where  A2  = A01  = Ai2  = A (X) , AiQ  = iA  U) , AQ0  = Q-A [X) , and 
A^  = Q-A (A+ip) , for  lii<c-l.  Aq  = cA  (tO  , and 
Aj_  = Q-A  (X+c\x)  . We  assume  that  jV>0_  and  u_>0_. 


Theorem  2 

Provided  p (c)  = (it A)  (cirp)  ^<1,  the  queue  is  stable. 
The  steady-state  vector  x = (xq  ,x1 , . . . ,2^,^  . . . ) is  given 

by 

(16)  X*  = xc_1r£"c+1,  for  k>c. 

The  matrix  Rc  is  the  unique  solution  in  the  set  of  non- 
negative matrices  of  spectral  radius  less  than  one  of  the 


! 

j 

; 

i 


J 

: 

i 
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equation 


(17)  cR^A(£)  + RC[Q-A  U+cy)  ] + A ( X)  =0. 

The  matrix  Rc  is  strictly  positive  and  cRcu.  = X. 
The  matrix  T,  given  by 


H 

O 

< 

o 

o 

10  An 

a12 

a20 

a21 

Ac-J 

[ Ac-1,0  Ac-l,l+  RcA0j 

is  an  irreducible  semi-stable  matrix  of  order  cm.  The 
vector  (xg  , x^ , . . . , x<=_^)  is  its  left  eigenvector  corresponding 
to  the  eigenvalue  zero.  It  is  normalized  so  that 

(3.9)  Xge  + ...  + xc_2e  + xc_1(I-Rc)  _1e  = 1. 

Proof 

The  assertions  about  the  matrix  Rc  were  already 
proved  as  part  of  Thm.  5 in  [1] . With  Rc  so  chosen,  the 
steady-state  equations 

— i-lA0  + *iAi  + xi+1A2  = 0,  i— c i 
are  satisfied,  provided  Xq,...,}^^  can  be  Pr°perly  chosen. 
The  c initial  equations  are  equivalent  to 
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Since  T is  clearly  irreducible,  it  suffices  to  verify  that 
T is  semi-stable.  The  off-diagonal  elements  of  T are  clearly 
nonnegative.  Moreover  we  have 

Ac-1 , 0—  + Ac-l,l-  + RcA0—  = RcA0—  “ Ac-1 , 2— 

= cRc_u  - X_  = 0^. 

It  follows  that  the  diagonal  elements  of  T are  negative  and 
hence  that  T is  semi-stable.  The  normalizing  condition  (19) 
uniquely  determines  x^,  0£i<c. 

Remark 

It  is  evident  that  the  proof  does  not  depend  on  the 
explicit  form  of  the  cmxcm  upper  left  hand  corner  of  Q* (c) , 
but  only  on  the  irreducibility  of  T.  Any  variations  of  the 
present  model  which  modify  only  the  upper  left  hand  corner 
will  therefore  have  a steady-state  vector  of  the  same  general 
! type. 

i 

m 
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IV.  A Numerical  Example  (c=l) 

In  order  to  illustrate  the  behavior  of  a queue  in 
which  short  periods  of  gross  instability  alternate  with 
long  periods  of  highly  undersaturated  conditions,  we  con- 
sider an  example  where  Q is  given  by  the  8x8  matrix 


-4 


4 0 0 0 


0 0 


0-1100 


0 0 


0 

0 


0 


Q = 


1 

4 


0 

0 

0 

0 

1 


0-110 


0 0 0 


00-110 


0 


0 


0 0 


0-1100 


9 


0 0 


00-110 


0 0 0 0 


0-11 


0 0 


0 0 0 


0 


-1 


and  the  vectors  X_  and  ^ are  given  by 


X = (40,  5,  2,  1,  1,  1,  1,1) , 

V = (10,  5,  5,  3,  3,  3,  3,  3)  . 

The  vector  tt  is  given  by 

. = £(1,  4,  4,  4,  4,  4,  4,  4), 
and  P = (irA)  (Try ) -1  = 0.8. 

*1  2 

The  equation  R = -A2A^  -R  AqA^  , was  solved  by 
successive  substitutions,  which  were  continued  until  the 
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maximum  entry-wise  difference  between  iterates  was  less  than 
— ft 

10  . This  took  989  iterations.  The  vector  was  evaluated 

and  it  was  found  that  < 10  ®e,  with  the  largest  differ- 

ence occurring  in  the  large  first  component. 


A small  improvement  in  accuracy  can  be  obtained  with 
little  additional  effort.  Let  R"±R'  be  the  last  two  iterates 
computed , then  we  compute  the  matrix  R by  Rj  ^ , =R^  ^ , + 

9 j (Rj j ' “Rj j • ) ' where  the  quantities  6j  are  chosen  so  that 
RA_  = 


The  conditional  means  and  variances  of  the  queue 


lengths  in  the  various  phases  were  computed  next  and  are 
listed  below. 


Table  1: 

Conditional 
the  Various 

Means  and  Varianc 
Phases 

Phase 

Mean 

Variance 

1 

50.65 

2575.55 

2 

50.98 

2581.12 

3 

41.55 

2482.10 

4 

36.05 

2346.10 

5 

31.29 

2178.66 

6 

27.17 

1996.99 

7 

23.60 

1812.20 

8 

20.50 

1631.56 

Selected  percentiles  of  the  conditional  distributions 
exhibit  very  clearly  the  variability  of  the  queue  length  over 
the  various  phases. 
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Table  2: 

Deciles 

of  the 

Conditional  Queue 

Length 

Distribution 

(the  first 

index  for 

which  : 

k%  is  exceeded) 

% 

Phase  1 

Phase  2 : 

Phase  3 

Phase 

4 Phase  5 Phase  6 

Phase  7 

Phase 

10 

5 

5 

0 

20 

11 

11 

2 

0 

0 

30 

18 

18 

8 

2 

1 

0 

0 

40 

26 

26 

16 

8 

2 

1 

1 

0 

50 

35 

35 

25 

18 

10 

4 

2 

1 

60 

46 

47 

36 

29 

22 

14 

7 

2 

70 

61 

61 

51 

44 

36 

29 

22 

14 

80 

82 

82 

71 

64 

57 

50 

42 

35 

90 

117 

117 

107 

99 

92 

85 

77 

70 

95 

152 

152 

142 

134 

127 

120 

113 

105 

99 

234 

234 

223 

216 

209 

202 

194 

185 

Table  3:  Conditional  Probabilities  of  Emptiness  in  the  Various 
Phases 


Phase 


Probability  of  Emptiness 


1 

2 

3 

4 

5 

6 

7 

8 


0.0142 

0.0168 

0.1284 

0.2079 

0.2701 

0.3233 

0.3692 

0.4090 


_ 
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Interpretation 

The  numerical  results  for  this  example  exhibit  a 
number  of  qualitative  features , which  deserve  to  be  stressed 
and  interpreted  as  they  are  not  present  in  more  elementary 
queueing  models. 

a.  The  global  mean  queue  length  is  given  by 
ir(I-R)  ^Re  = 33.63,  but  this  is  clearly  not  a 
meaningful  descriptor  of  this  highly  oscilla- 
tory queue. 

88 

b.  The  M/M/1  queue  with  parameters  A°=TrX=Tg-,  and 

110  t 

la  -Tnj.=— » has  a stationary  mean  queue  length 
of  4,  with  a variance  of  20.  It  does  not  begin 
to  offer  an  approximation  to  the  present  queue. 

c.  The  parameters  chosen  for  the  numerical  example 
can  be  thought  of  as  representing  a short  over- 
saturated rush  hour  (Phase  1) , a transitional 
period  (Phases  2 and  3)  and  an  unsaturated  period 
(Phases  4 through  8) , repeated  cyclically.  The 
conditional  queue  length  distributions  at  the 
ends  of  Phases  1,3  and  8 represent  the  most 
"extreme"  queue  conditions.  By  mixing  the  con- 
ditional queue  length  densities  of  the  Phases  4 , 
5,..., 8 with  the  weights  it ^ (ir^+. . .+7Tg)  \ 4£j£8, 
we  obtain  the  conditional  queue  length  density, 
given  that  the  queue  is  in  its  unsaturated  period. 
This  density  will  clearly  be  different  from  that 
at  the  end  of  Phase  8. 
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d.  Although  the  mean  queue  lengths  vary  consider- 
ably between  phases  1 and  8,  there  is  much  less 
variation  in  the  standard  deviations.  This  is 
due  to  two  causes.  During  the  oversaturated 
phase  there  is  enormous  random  variability  in  the 
behavior  of  the  queue.  Very  long  queues  at  the 
end  of  a rush  hour  will  take  a long  time  to  dis- 
sipate, while  shorter  queues  dissipate  quickly. 
This  effect  is  likely  to  become  apparent  in  high 
variances  in  spite  of  reduced  means  in  the  higher 
phases. 

The  second  cause  of  the  high  variability  in  those 
phases  lies  in  the  exponential  distribution  of 
the  sojourn  times  in  each  phase.  The  effect  of 
this  assumption  is  easy  to  study.  By  modifying 
the  Q-matrix,  we  can  change  the  distributions  of 
the  rush  hour,  the  transition  period  and  the 
undersaturated  period  in  a versatile  manner  with- 
out changing  their  mean  durations.  These  periods 
can  have  arbitrary  phase  type  distributions, 
which  may  have  arbitrary  positive  coefficients  of 
variation. 

In  order  to  illustrate  the  latter  points  numerically  we  com- 
puted a number  of  parametric  variants  of  this  model.  The 
results  will  be  presented  only  briefly,  but  they  are  striking 
indeed . 
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Variant  1 

We  replace  the  exponential  duration  of  the  rush-hour 
by  an  Erlang  distribution  of  order  2 with  parameter  2.  This 
does  not  affect  the  mean  duration  of  the  rush-hour  period, 
but  reduces  its  variance.  The  inner  products  ji_x  and  try 
remain  as  before.  There  are  now  nine  phases. 


Table  4:  Conditional  Means  and  Variances  in  the 
Various  Phases 


Phase 

Mean 

Variance 

0 

26.69 

801.99 

1 

41.70 

1051.32 

2 

41.87 

1077.30 

3 

32.03 

1061.05 

4 

26.34 

985.76 

5 

21.53 

883.98 

6 

17.53 

771.76 

7 

14.24 

660.04 

8 

11.54 

555.57 

Phases  0 and  1 now  correspond  to  the  rush-hour.  We  note  the 
drastic  reduction  in  the  means  and  variances  in  all  the 
phases.  The  dependence  on  the  random  variability  of  the 
duration  of  the  rush-hour  period  is  very  strong  indeed. 

Variant  2 

In  this  case  we  increased  the  service  rate  during  the 
transitional  period,  keeping  iry  constant  by  reducing  the  ser- 
vice rate  later  in  the  undersaturated  case.  The  parameters 
are  as  in  the  original  model,  except  for  y.  which  is  now 
(10,  7,  7,  3,  3,  3,  1,  1).  The  effect  of  this  change  in 
parameters  is  not  as  substantial  as  one  might  expect. 
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Table  5 : Conditional  Means  and  Variances  in 
the  Various  Phases 


Phase 

Mean 

Variance 

1 

54.56 

2775.76 

2 

47.94 

2734.02 

3 

34.95 

2442.26 

4 

30.49 

2258.57 

5 

26.60 

2067.01 

6 

23.21 

1875.74 

7 

23.98 

1846.81 

8 

24.51 

1828.46 

Variant  3 

The  most  striking  improvement  in  the  conditions  of 
the  queue  is  obtained  by  increasing  the  service  rate  during 
the  rush  period  itself.  This  is,  of  course,  not  always 
economically  feasible. 

The  parameters  are  again  as  in  the  original  model , 
except  for  v which  is  now  (30,  5,  5,  2,  2,  2,  2,  2). 


Table  6:  Conditional  Means  and  Variances  in 
the  Various  Phases 


Phase 

Mean 

Variance 

1 

15.75 

251.15 

2 

16.66 

260.63 

3 

10.13 

214.98 

4 

8.49 

187.23 

5 

7.12 

161.06 

6 

5.99 

137.07 

7 

5.07 

115.67 

8 

4.31 

97.01 

A final  comment  is  concerned  with  the  relationship 
between  the  means  and  corresponding  standard  deviations  in 
all  the  preceding  examples.  We  see  that  queues,  which 


t 


_ . 
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periodically  pass  through  phases  in  which  they  are  locally 
oversaturated,  exhibit  considerable  random  variation.  To 
obtain  reliable  parameter  estimates  by  Monte  Carlo  simula- 
tion, with  the  rather  small  sample  sizes  employed  in  prac- 
tice, appears  to  be  impossible.  The  numerical  behavior  of 
the  present  examples  casts  serious  doubt  on  the  merits  of 
simulation  methods,  when  applied  to  realistic  queueing 
models,  unless  prohibitively  large  sample  sizes  are  used. 
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